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EXISTENCE OF PERIODIC SOLUTIONS OF THE 
FITZHUGH-NAGUMO EQUATIONS FOR AN EXPLIGIT RANGE 
OF THE SMALL PARAMETER* 

ALEKSANDER CZECHOWSKit AND PIOTR ZGLICZYNSKit 


Abstract. The FitzHugh-Nagumo model describing propagation of nerve impulses in axon is 
given by fast-slow reaction-diffusion equations, with dependence on a parameter e representing 
the ratio of time scales. It is well known that for all sufficiently small e > 0 the system possesses 
a periodic traveling wave. With aid of computer-assisted rigorous computations, we prove the 
existence of this periodic orbit in the traveling wave equation for an explicit range e E (0, 0.0015]. 
Our approach is based on a novel method of combination of topological techniques of covering 
relations and isolating segments, for which we provide a self-contained theory. We show that 
the range of existence is wide enough, so the upper bound can be reached by standard validated 
continuation procedures. In particular, for the range e E [1.5 x 10”'^, 0.0015] we perform a rigorous 
continuation based on covering relations and not specifically tailored to the fast-slow setting. 
Moreover, we confirm that for e = 0.0015 the classical interval Newton-Moore method applied to a 
sequence of Poincare maps already succeeds. Techniques described in this paper can be adapted to 
other fast-slow systems of similar structure. 

Key words, fast-slow system, periodic orbits, rigorous numerics, FitzHugh-Nagumo model, 
isolating segments, covering relations 
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1. Introduction 


( 1 ) 


1.1. The FitzHugh-Nagumo equations. The FitzHugh-Nagumo model with 
diffusion 

du 1 d'^u , ... . 

— = + u(u - a)(l -u) -w, 

OT 7 ox^ 

dw 

— =e{u-w), 

was introduced as a simplification of tfie Hodgkin-Huxley model for the nerve impulse 
propagation in nerve axons 13 The variable u represents the axon membrane 
potential and w a slow negative feedback. Traveling wave solutions of (Hare of 
particular interest as they resemble an actual motion of the nerve impulse |20| . By 
plugging the traveling wave ansatz (u,w){t,x) = {u,w){x + Or) = {u,w){t), 0 > 0 
and rewriting the system as a set of first order equations we arrive at an ODE 


u = V, 

v' = ^{9v — u{u — a)(l — It) -I- w), 


( 2 ) 


w' = -{u — w). 

U 


to which we will refer to as the FitzHugh-Nagumo system or the FitzHugh-Nagumo 
equations. The FitzHugh-Nagumo system is a fast-slow system with two fast 
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institute of Computer Science and Computational Mathematics, Jagiellonian University, 
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(a) Approximate periodic orbit for (^) Projection onto (u, w) plane. 

£ = 0.001 in blue. 
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(c) A schematic drawing of the singular 
periodic orbit. 

Figure 1. A numerical approximation of the periodic orbit close 
to the singular orbit, the slow manifold in green. 


variables u^v and one slow variable w. The parameter 9 represents the wave speed 
and e is the small parameter, so 0 < e <C 1. To focus our attention, following 4|[6p8 
we set the two remaining parameters to 


a := 0.1, 


7:= 0.2, 


( 3 ) 


throughout the rest of the paper. 

Bounded solutions of ([^ yielding different wave profiles have been studied by many 
authors both rigorously and numerically, see [3l[4}[t|[^[8| [T^[T^[l^[^[^[M}|^|30[|^ 
and references given there. Periodic orbits leading to periodic wave trains exist 
for an open range of 0’s and were treated in [^[^[^[^[^[^ ; traveling pulses 
generated by homoclinic orbits exist for two isolated values of 0, their existence was 
proved in [^|4l[^[^ . Stability of waves was discussed in . Proofs of 


existence use various methods, but most share the same perturbative them^ We 
outline it below, for the periodic orbit. 

Consider the limit equation at e = 0. There, the velocity of w is zero and 
the phase space can be fibrated into a family of two-dimensional fast subsystems 
parameterized by v. These subsystems serve as a good approximation to the system 


bn 


3 [4 the authors use non-perturbative computer-assisted methods for a single value e = 0.01 


where the system becomes a regular, although stiff ODE. It is to be noted that in authors prove 
stability of the wave for this particular parameter. 
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with e > 0 small, except for regions of phase space near the Z-shaped slow manifold 


Co = {(w,n,u>) : v = 0, w = u(u — 0.1)(1 — m)}, 


(4) 


where the velocities of fast variables become small and the slow variable takes 
over. There, one can analyze the slow flow given by evolving w while keeping the 
restriction (m, n) € Co, obtaining a differential algebraic equation. Note that the 
slow manifold is formed by fixed points of the fast subsystems. For a range of w 
it has exactly three branches - by looking from a perspective of the respective u 
values the lower and the upper one are formed by saddles, and the middle one is 
formed by sources. We denote the upper/lower branches of saddles by A„(w) and 
Ad{w), respectively. For exactly two values w G {w^,,w*}, with w* < w* there are 
heteroclinic connections from Ad(rt;,) to A„(u;*) and from A„(u;*) to Aflw*). To 
see that it is best to regard the fast subsystem as a Hamiltonian system with added 
negative friction and shoot with w, a detailed description is given in [^. It happens 
that in the range [w*, w*] the slow flow on the branch A^^ is monotonically decreasing 
and on the branch A^ monotonically increasing, so by connecting the heteroclinics 
with pieces of the slow manifold one assembles the singular periodic orbit, see 
Figure The proof of existence of an actual periodic orbit goes by perturbing to 
e > 0 small and using certain arguments based on topological methods or 


Fenichel theory and differential forms 23,24 


Using methods described above proofs of existence have been given for e G (0, Eq]) 
Eo “small enough”. With aid of computer we are able to improve these results. We 
give an explicit eq such that a periodic orbit of ([^ exists for the parameter range 
E G (0, Eo]. Our secondary objective is to make eq as large as possible, so there is no 
doubt that for e > eq one can perform further continuation using well-established 
computer-aided methods such as the interval Newton-Moore method [2 29,31 
applied to a sequence of Poincare maps. The main results of this paper are: 


Theorem 1.1. For each e G (0,1.5 x 10 ^], for 9 = 0.61 and other parameter values 
given in ([^ there exists a periodie orbit of ([^ . 

Theorem 1.2. For each e G [1.5 x 10“^, 0.0015], for 6 = 0.61 and other parameter 
values given in ([^ there exists a periodie orbit of ([^. 

Theorem 1.3. For e = 0.0015, 9 = 0.61 and other parameter values given in ([^ 
there exists a periodic orbit of <I1|. which is formed from a locally unique fixed point 
of a Poincare map. 


The reason we do not merge statements of Theorems 1 1.1 1 and 1.2 is a significant 
difference in proof techniques. For the proof of Theorem |1.1| we exploit the fast-slow 
structure and construct a sequence of isolating segments and covering relations 
around the singular orbit. Let us remark that combination of these two topological 
methods is new, and isolating segments have been previously applied only to time- 
dependent equations. For Theorem 1.2 we perform a “regular ODE” type of proof, 
with a parameter continuation method based on verifying covering relations around 
an approximation of the periodic orbit. Finally, in Theorem |1.3| we are far enough 
from E = 0, so that a proof by the interval Newton-Moore method applied to a 
sequence of Poincare maps succeeds (see Subsection 4.4), establishing both the 
existence and local uniqueness. 

The motivation for this choice of wave speed was that in numerical simulations 
parts of the periodic orbit near the slow manifold stretched relatively long, which 
allowed us to fully exploit its hyperbolicity. In the program files a lot of values 
were hardcoded for this particular 9, but we report that by substituting 9 = 0.53, 
9 — 0.47, 9 G [0.55,0.554] we were also able to produce results like Theorem [m for 
a (shorter) range of e G (0, 5 x 10“®]. We think that by spending time tuning the 
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values in the proof, the range of e for these 0’s could have been made wider. This is 
of course the very same orbit and if one had enough patience, then continuation in 
0 would be theoretically possible. 

Let us observe that proofs employing interval arithmetics are fairly easy to 
adapt to compact parameter ranges. A validated continuation can be performed 
by subdividing the parameter interval finely and feeding the program with small 
parameter intervals instead of one exact value 11 14 26 36 . However, we emphasise 


we are dealing with a half-open 


again that when trying to prove Theorems |1.1 
range and such straightforward approach is bound to fail - as e —>■ 0+ the period 
of the periodic orbit grows to inhnity and in the singular limit e = 0 the orbit is 
destroyed. Therefore not all premises in theorems implying existence of such orbits 
can be verified for e G [0, eo]; the assumptions need to be formulated in such a 
way, that the ones which are computationally difficult are possible to check with 
computer aid for e € [0, Cq]: ^-nd the leftover ones are in a simple form where e can 
be factored out “by hand”, assuming e > 0. To be more specific, e is factored out to 
check the monotonicity of the slow flow in certain sets (the aforementioned isolating 
segments) in the neighborhood of the slow manifold, which reduces to a simple 


inequality for the slow part of the vector field (inequality (Sla) in Remark 2.16 cf. 


inequality (34)). Clearly, such technique can be applied to other fast-slow systems 
with one-dimensional slow manifold, even if the slow flow is nonlinear. 

Full proofs are executed with computer assistance and described in detail in 
Section!^ For the computer assisted assumption verification we use the CAPD 
library [n, which provides algorithms for enclosures of solutions of ODEs, Poincare 
maps and their derivatives. By a rigorous enclosure we mean that the true result of 
a given operation is always contained in an interval object (vector, matrix) returned 
by the program procedure. Since our assumptions are in essence a collection of 
strict inequalities, if we keep the overestimates small, we should be able to verify 
them on a machine. 

The code which executes the necessary computations is available at the author’s 
homepage [^. Below we outline the basic ideas of each proof. 

1.2. Outline of the proof of Theorem |1.1[ We conduct a phase space proof 
based on a reduction to a sequence of Poincare sections and a fixed point argument 
for a sequence of Poincare maps. For e > 0 small the orbit switches between two 
regimes - the fast one close to heteroclinics of the fast subsystem and the slow one 
along the branches A^, A„ of the slow manifold. The strategy is to form a closed 
sequence of covering relations and isolating segments and deduce the existence of a 
fixed point of a sequence of Poincare maps via a topological theorem - Theorem |2.30[ 

For the fast regime we employ the rigorous integration and we check covering 


relations as defined in 39 . Informally speaking, one introduces compact sets called 
h-sets on the Poincare sections near the points h-d{w<,), h.u{w^), Au{w*), Ad{w*). 
They come equipped with a coordinate system with one direction specified as exit, 
the other as entry. To verify a covering relation by a Poincare map between such two 
sets X, Y one needs to check that the exit direction edges of X are stretched over 
Y in the exit direction and that the image of X is contained in the entry direction 
width of Y, see Figure in Section The “shooting” in the exit direction is in fact 
made possible by a non-degenerate intersection of stable and unstable manifolds of 
the respective fixed points in the singular limit e = 0, as described in condition (P2) 
in Section HI 

Around the slow manifold branches we place isolating segments, which allow us 
to track the orbit in this region without rigorous integration, see Theorem |2.22| In 
short, isolating segments are three-dimensional solids diffeomorphic to cubes and 
akin to isolating blocks of Conley and Easton or periodic isolating segments of 
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Srzednicki 32p3 . For each isolating segment one distinguishes three directions: exit, 
entry and central. It is required that the faces in the exit direction are immediate 
exit sets for the flow, the faces in the entry direction are immediate entrance sets 
and the flow is monotone along the one-dimensional central direction. The first 
two assumptions are checked by a computer (for e S [0,1.5 x 10“^]), exploiting the 
hyperbolicity of branches of the slow manifold. The last one we can easily fulfill 
by aligning the central direction of segments with the slow variable direction. This 
setup reduces the central direction condition to a verification of whether ^ 0 for 

all points in such segment. Under assumption e 0 we can then factor out e from 
the slow velocities, and our condition reduces to a question whether u ^ w for all 
points in each segment, which is straightforward to check. This is the only moment 
in the proof, when we need to make use of the assumption that e is strictly greater 
than 0 (i.e. e S (0,1.5 x 10“^]). 

We place four supplementary “corner segments”, containing the corner points 
Kd{wt), A„(in*), Au{w*), Ad{w*), the role of which is to connect the h-sets with the 
segments around the slow manifold. From the viewpoint of definition these are no 
different than regular isolating segments. However, the mechanism of topological 
tracking of the periodic orbit here is slightly distinct, as the central direction changes 


roles with the exit/entry ones, see Theorems 2.28 2.29 


To obtain a closed loop, the sizes of the first and the last h-set in the sequence 
have to match. For that purpose isolating segments around the slow manifold may 
need to grow in the exit direction and compress in the entry one as we move along 
the orbit. This way we can offset the size adjustments of the h-sets, which may be 
necessary to obtain covering relations in the fast regime. The analysis of a model 
example performed in Theorem |3.1| is devoted to providing an argument for why 
this should work for e small. The main idea is that as e —> 0'^ the vector field in 
the slow/central direction decreases to 0 and the slope of the segment becomes 
irrelevant when checking isolation, see Figure in Section 

A schematic drawing representing the idea of the proof for the model example is 
given in Figure]^ in Section]^ 


1.3. Outlines of the proofs of Theorems 1.2, 1.3 We are already at some 


distance from e = 0, but for small e the periodic orbit’s normal bundle is consisting 
of one strongly repelling and one strongly contracting direction, so any attempts of 
approximating the orbit by numerical integration, either forward or backward in 
time, fail. On the other hand, the singular orbit at e = 0 no longer serves as a good 
approximation for the purpose of a computer assisted proof. As we can see, the 
challenge now is on the numerical, rather than the conceptual side. To find our good 
numerical guess we introduce a large amount of sections, so that the integration 
times between each two of them do not exceed some given bound and then apply 
Newton’s method to a problem of the form 


Pi{xi) -X 2 = 0, 

P2{x2) - a;3 = 0, 

Pk{xk) -xi =0, 


(5) 


where Pi’s are the respective Poincare maps. Then, we construct a closed sequence 
of h-sets on these sections and verify covering relations between each two consecutive 


ones, to prove the periodic orbit by means of Theorem 2.8 (Corollary 7 in 39 ). 
Since we control the integration times, isolating segments are not needed anymore - 
in Theorem o they were used for pieces of the orbit where the integration time 
tended to infinity. 
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By a rigorous continuation with parameter e, we are able to get an increase 
of one order of magnitude for the upper bound of the range of e’s, for which the 
periodic orbit is confirmed. Without much effort we show that for this value of e 
the classical (see and references given there) method of application 

of the interval Newton-Moore operator to a problem of form ([^ succeeds. This 
requires a rigorous computation, but these are handled efficiently by the 


Lohner algorithm implemented in CAPD 38 


1.4. Organization of the paper. The content of this paper is arranged as follows. 
In Section we provide a self-contained theory on how to incorporate isolating 
segments into the method of covering relations. Typical application to find periodic 


orbits is given by Theorem 2.23 which works for any number of exit and entry 
dimensions. This theorem is however only for future references; for the FitzHugh- 
Nagumo equation we use Theorem |2.30[ which is restricted to one exit and one 
entry direction and employs a certain switch between slow and fast directions in 
segments (Theorems |2.28 |2.29 ), which can be viewed as a topological version of 
the Exchange Lemma (cf. Chapter 5 in [23] and references given therein) from 
geometric singular perturbation theory (GSPT). In Section|^we give arguments of 
why one can find a closed sequence of covering relations in the FitzHugh-Nagumo 
system for e > 0 small: given a 3D model fast-slow system with a singular orbit 
sharing the qualitative properties of the one in ([^, we argue that it is possible to 
apply Theorem |2. 30 to obtain the existence of a periodic orbit. In some sense this 
is a repetition of methods in (^[8 20 recast in the language of covering relations 


and isolating segments. The exposition of Section is informal and no part of this 
section is needed to prove any theorems outside it, in particular Theorems 1 1.1 [ [L^ or 
1.3| the reason for including it in this paper is to give the reader some intuition for 
the relation between our abstract topological methods and the fast-slow structure 
of the system. Finally in Section]^ we give details of the computer assisted proofs 
of Theorems |1.1[[L^|1.3| and provide some numerical data from the programs. 


2. Covering relations and isolating segments 

2.1. Notation. Unless otherwise stated H-H can be any fixed norm in K". Given a 
norm, by i?„(c, r) we will denote the ball of radius r centered at c G M". We will 
drop the subscript if the dimension is clear from the context. By (•, •) we will denote 
the standard dot product in K". 

We assume that K is always equipped with the following norm: ||a;|| = \x\. 
Given a set Z, by int Z, Z, dZ and conv Z we will denote the interior, closure, 
boundary and the convex hull of Z, respectively. 

Given a topological space X, a subspace D c M. x X and a local flow ip : D ^ X, 
by writing ip{t,x) we will implicitly state that {t,x) G D, so for example by 

(p{t,x) = y (6) 

we will mean (/?(t, x) exists and (p{t, x) = y. 

By idjc we denote the identity map on X. The symbol const denotes a constant - 
usually some uniform bound - the value of which being not important to us, so for 
example the expression f{x) > const, x G X, means that there exists C G M such 
that f{x) > C yx G X. 

The local Brouwer degree of a continuous map g : K" -G- M" at c G M" in an open, 
bounded set fl C R" will be denoted by deg{c,g, fl), provided it is well-defined. 

By smoothness we mean smoothness. In some assumptions differentiability 
would be enough, but we do not go into such details. 
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2.2. Covering relations - basic notions. In this subsection we recall the defini¬ 
tions of h-sets, covering and back-covering relations for maps as introduced in 


39 


We make a following change in the nomenclature: in 39 various objects related to 
h-sets (directions, subsets, etc.) are being referred to as unstable or stable. We will 
refer to them as exit and entry/entrance, respectively; we think that this reflects 
better their dynamical nature and does not lead to misunderstandings. However, we 
keep the original symbols u, s, so u should be connoted with exit and s with entry. 


Definition 2.1 (Definition 1 in [^). An h-set is formed by a quadruple 

X = {\X\,u{X),s{X),cx) (7) 

consisting of a compact set \X\ C K" - the support, a pair of numbers u(A'), s(A) G N 
such that u{X) + s(A) = n (the number of exit and entry directions, respectively) 
and a coordinate change homeomorphism cx '■ K” —>■ x such that 

cv(|X|) = H„(;f)(0,l)xH,(x)(0,l). (8) 

We set: _ _ 

Xc ■= -Bm(x)(0! 1) X 7?s(a)(0: 1)5 

:= dB^{x)iO. 1) X i?,(x)(0,l), 

X+ := H„(x)(0,l) X 9H,(x)(0,1), (9) 

A- :=c^i(A-), 

A+ :=cj^i(A+). 

We will refer to X~/X^ as the exit/entrance sets, respectively. To shorten the 
notation we will sometimes drop the bars in the symbol |A| and just write X to 
denote both the h-set and its support. 


Remark 2.2. Due to condition ([^, it is enough to specify u{X),s{X) and cx to 
unambiguously define an h-set X. 

Definition 2.3 (Definitions 2, 6 in [^). Assume that X,Y C M" are h-sets, such 
that u{X) = u{Y) = u and s{X) = s(T) = s. Let 5 : D —)■ K" be a map with 
|A| C D C M". Let gc = cy o g o : X^ —>■ IR“ x M'* and let w be a non-zero 
integer. We say that X g-covers Y with degree w and write 

A ^ A (10) 

iff g is continuous and the following conditions hold 

1. there exists a continuous homotopy h : [ 0 ,1] x A^ —>■ x M®, such that 

ho=gc, ( 11 ) 

/i([O,l],A-)nFe = 0 , (12) 

M[0,l],A,)nr+ = 0 . (13) 

2. There exists a continuous map A : ]R“ —>■ ]R“, such that 

hi{p,q) = (A(p),0) Vp G H„(0,1) and q G 5^(0,!), 
A(aH„(0,l))cK“\i?„(0,l), (14) 

deg( 0 , A,H„( 0 , 1 )) = w. 

In case A is a linear map, from (A 8 ) we get deg (0, A, i?„(0,1)) = sgndet A = ± 1 . 
In such situation we will often say that A p-covers Y, omit the degree and write 

A=^r. 


Remark 2.4 (Remark 3 in [^). For m = 0 we have 5^(0,1) = 0 and A p-covers Y 
iff g(|A|) is a subset of int |A|. In that case, we formally set the degree w to 1. 
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Definition 2.5 (Definition 3 in 39 ). Let X be an h-set. We define the transposed 
h-set as follows: 

• u{X'^) = s{X) and s{X^) = u{X), 

• Cxt{x) = j{cx{x)), where j : x —)• x is given by 

jip,q) = iq,p)- 

Observe that {X’^)~^ = X~ and {X'^)~ = X'^, thus transposition changes the 
roles of exit and entry directions. 


Definition 2.6 (Definition 4, 7 in 39 ). Let X,Y be h-sets with u{X) = u{Y) and 
s(X) = s{Y). Let g : fl C M" —t K". We say that X ^-backcovers Y with degree w 
and write X Y iff 

• g-^ ■■ 

• g- 


1 ” exists and is continuous, 
'’“^-covers with degree w. 


Definition 2.7 (Definition 5 in 


39 


We will use the notation X 


.g,™ 


Y and say 


that X generically ^-covers Y with degree w iff any of these two hold: 

• X (/-covers Y with degree w, 

• X (/-backcovers Y with degree w. 

Again, we will sometimes omit the degree in our notation, in case the homotopy 
can be given to a linear map. 

Let us state the fundamental lemma motivating the use of covering relations for 
finding trajectories and fixed points of sequences of maps. 

Theorem 2.8 (Theorem 9 in [^). Let Xi, f S {0,..., /c} be h-sets with u{Xq) = 
• • • = u{Xj.), s(Aro) = • • • = s{Xk) and set n = u{Xq) -|- s(Xo). Assume that we 
have the following chain of covering relations: 

gk,Wk, 


TA AA .92,'^2. -VA .93,'^3. incite, AA 

Aq < ■ ■ ■ > Ai <■■■ ■ > A2 <■■ - > . . . < ■ - ■ > Afe. 
for some Wi S Z*. Then there exists a point x G int Xq such that 
(Piogi-i o ■ ■ ■ogi){x) e intx^, f e {i,... k}. 
Moreover, if Xk = Xq, then x can be chosen so that 

{gk o gk-i o • • • °gi){x) = x. 

For an h-set X with u(X) = 1, s(X) = s we have: 


W = [-l,l]xB,(0,l), 

= ({-1} X u ({1} X 


(15) 

(16) 

(17) 

(18) 


We will often use the following geometrical criterion for verifying of the covering 
relation in such a case: 


Lemma 2.9 (Theorem 16 in 39 ). Let X,Y be h-sets with u{X) = u(Y) = 1, 
s(X) = s(Y) = s. Let g : X ^ be a continuous map. Assume that both of the 

following conditions hold: 

(Cl) We have 

gc{Xc) C int(((—oo, —1) x M®) U Yc U ((1, oo) x K®)), (19) 

(C2) either 


5c({-l} X i?s(0, 1)) C (-00, -1) X K® and gc({l} x i?s(0, 1)) C (1, oo) x 


(/c({ —1} X Bs{0, 1)) C (1, oo) X K® and gc({l} x Bs{0, 1)) C (—oo, —1) x M®. 

( 20 ) 
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{CY o g){X 



Figure 2. A covering relation X Y. 


Then 

X =M=> Y. (21) 

Remark 2.10. In applications it is convenient to introduce the notation X~’^ = 
1} X Bs{0, 1)) {the left exit edge) and X~’^ = c)^^({l} x i?s(0, 1)) {the right 
exit edge) and check (Cl), (C2) by putting 

9c{Xc) = {cYog){\X\), 

g,{{-l}xB;{^) = {cYog){X-’^), (22) 

5 c({l}x(B;(^) = (cy0 5)(A-’’'), 

see Figure 

Analogously, if for an h-set Y we have u{Y) = u and s(F) = 1, then 
re = :B:(Mx [-1,1], 

y/ = (il„(0,l) X {-!}) U (i?„(0,i) X {!}), 
and we can apply the same principle to transposed sets: 

Lemma 2.11. Let X,Y be h-sets with u{X) = u(Y) = u and s(A) = s(y) = 1. 
Let LI C ]R“+^ and g : LI ^ K“+^ be continuous. Assume, that g~^ : |y| —>■ 
exists, is continuous and that both of the following conditions hold: 

(Cl a) We have 

gf^iYc) C int((M“ x (—oo, —1)) U Ac U (IR“ x (1, oo))); (24) 

(C2a) either 

gf^{Bu{0, 1) X {-1}) C K” X (-00,-1) and gf^{Bu{0, 1) x {!}) C K“ x (l,oo) 

or 

gf^{Bu{0, 1) X {—1}) C X (1, oo) and gf^{Bu{0, 1) x {!}) C K” x (—oo, —1). 

(25) 

Then 

X Y. (26) 

In such case we will sometimes operate with the notation Y^^ = Cy^(i?„(0,1) x 
{ —1}) {the left entrance edge) and F+t = Cy^(i?u(0,1) x {1}) {the right entrance 
edge). Conditions (Cla) and (C2a) can then be rephrased in the same manner as in 
Remark 12.101 
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2.3. Covering relations and Poincare maps. In this subsection we will consider 
an ODE 


i = f(x), 


( 27 ) 


given by a smooth vector field /, and describe how h-sets and covering relations 
can be used in such setting. 

Assume, that we are given a diffeomorphism $ : —?> , and let E C be 

a subset of the hypersurface S ;= $“^({0} x A point a; S S is regular iff 

{n{x), f{x)) ^ 0, where n{x) is a normal to S at x. If every point a; G S is regular, 
then we will say that E is a transversal section. 

For a given xq G we will denote by 


Vj(<,a:o) 

the local flow generated by /, that is the value of the solution 
initial condition a;(0) = Xq. 

Let El, E 2 be two transversal sections such that: 


(28) 

x(t) to (27) with the 


• either Ei C E 2 or Ei n E 2 = 0, 

• we have int Ei = Ei and int E 2 = E 2 , where closures and interiors are taken 
in the hypersurface topology, 

• for each x G Ei there exists a r > 0 such that (/?(t, x) G E 2 
It is well known that the Poincare map: 


P : El 9 X 


inf <p(t, x) G E 2 
T:<^(r,£c)GS2 


(29) 


is well-defined and smooth for points x G int Ei such that P(x) G int E 2 (interiors 
in the hypersurface topology). The proof can be found in e.g. [^ . 

To make the formulation of some theorems in future shorter, we extend the above 
definition of a Poincare map in the scenario Ei C E 2 to also cover the embedding 
by identity id : Ei —>■ E 2 . In such case we will always specifically refer to such map 
as the identity map, to differentiate from a Poincare map P given by (29). 

For such a Poincare map we define the h-sets in a natural manner. We can 
identify Ei, E 2 with two copies of Then we can proceed to describe the h-sets 

on each of these copies - note that they will be h-sets in not K^. 


Remark 2.12 . Tr eating h-sets as subsets of sections is a slight abuse when compared 
to Definition 2.1 where they were subsets of the Euclidean space Nevertheless, 


we can always compose the change of coordinates homeomorphism for the h-set 
with the global coordinate frame on the section to get back to the Euclidean space. 
Therefore, given a section E C for an h-set A C E the actual coordinate 
change will take the form cx = cx o where $ : E —>■ {0} x is the global 

coordinate frame for the section and cx : —?> x is a coordinate 

change homeomorphism satisfying ([^. 

2.4. Isolating segments. Assume we are given a smooth vector field (27), an 
associated local flow (28) and a pair of transversal sections Ej^, Eout- 


Definition 2.13. A segment between two transversal sections Ejn and Eout is formed 
by a quadruple S = (l^l, u(S'), s(S'), eg), consisting of a compact set [S'] C {the 
support), a pair of numbers u{S), s{S) G N with u{S) + s{S) = N —1 (the number of 
exit and entrance directions, respectively) and a coordinate change diffeomorphism 
eg : —>■ X X M such that: 


cg(| 5 |) =i 3 „(g)( 0 ,l) X i 3 ,(g)( 0 ,l) X [ 0 , 1 ], 
Cg^(^«(S)( 0 j 1 ) X -Ss(g)(0, 1 ) X { 0 }) C Ein, ( 30 ) 

Cg^(^«(S)(0i 1) X i3s(g)(0, 1) X {!}) c Eout- 
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We set: 

Sc ■= S„(s)(0,1) X 1) X [0,1], 

5- := 35,(5) (0,1) X 5,(5) (0,1) x [0,1], 

5+ := 5,(5) (0,1) X 35,(5) (0,1) x [0,1], (31) 

5 - ■.= c-s\S:), 

S^ -=€3^8+). 

We will refer to S~/S~^ as the exit/entrance sets, respectively. Again, to shorten 
the notation sometimes we will drop the bars in the symbol |5| and just write S to 
denote both the segment and its support. 

Remark 2.14. As with h-sets, it is enough to give u{S), s{S) and 0$ to define a 
segment S. 

Given a segment S we introduce the following notation for projections: 


7r„ : X x M 9 (x,,a;,,a;^) —>■ x, G 

TT, : X X K 9 {xu, a;,, a;^) ^ x, G (32) 

TT^ : X X M 9 {xu,Xg,x^) — >■ G M. 


Definition 2.15. We say that S is an isolating segment between two transversal 
sections Ej, and Eout if 5 is a segment, the functions x —>■ ||7r„(x)||, x —>■ ||7r,(x)||, 
X ||7r^(x)||, X G X X M are smooth everywhere except at 0 and the 

following conditions are satisfied: 

(51) )^7r^C5(i)3(f, x))|j^p > 0 for all x G [S'! {monotonicity), 

(52) ||7r„C5(:p(t, x))|||^_^ > 0 for all x G S~ {exit set isolation), 

(53) ||7rsC5(:p(t,x))|||^_^ < 0 for all x G S~^ {entrance set isolation). 


As one can see, our definition of an isolating segment S relies on splitting the 
phase space into: 


• the exit directions 71 , 005 , 

• the entry directions 71 , 005 , 

• the one-dimensional central direction tt^ o 05. 


In that sense, it is a simplification of the concept of periodic isolating segments 

(also, under a name 


in nonautonomous systems, as originally introduced in 
of periodic isolating blocks in [32| ), where a wider range of boundary behavior 
was considered. On the other hand, contrary to 33 , we are able to work with an 


autonomous ODE - in 33 the central direction had to be given by time. 

When introducing an isolating segment we will sometimes omit specifying the 
transversal sections, as they are implicitly defined by C5. 


Remark 2.16. Each of the conditions (S1)-(S3) is equivalent to its following counter¬ 
part: 

(Sla) (V( 7 r^ ocs){x),f{x)) > 0 for all x G 151, 

(S2a) (V IItt, o C 5 II (x), /(x)) > 0 for all x G S~, 

(S3a) (V IItt, o 05 ] | (x), /(x)) < 0 for all x G 5+. 

Since S~, S~^ are subsets of the level sets {x G : ||7r„ o C 5 (x)|| = 1}, {x G 
: Iks o C 5 (x)|| = 1}, respectively, and gradients are normals to level sets, (S2a) 
and (S3a) can also be rewritten as: 

(S2b) {n-{x), f{x)) > const > 0 for all x G S~, 

(S3b) k+(x), /(x)) < const < 0 for all x G 5+, 
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respectively, where n::p(x) are normals to pointing in the outward direction of 

In our applications the faces of segments will always lie in affine subspaces, hence 
conditions (S2b) and (S3b) are easy to check by an explicit computation. Let be 
the projection onto i-th variable, i G 1.,..., N. In the central direction our changes 
of coordinates will take an affine form 

T^fj.cs{x) = aTTi{x) + b, (33) 

for a 7 ^ 0, b G K. In that situation (Sla) is equivalent with 

sgn(a)7r,(f(x)) > 0, Vx € l^l, (34) 

which again is easily established. In particular, if the sign of 7ri(f(x)) is negative, one 
needs to orient the segment in the direction reverse to the z-th coordinate direction 
by giving a a negative sign. 

We will now introduce the notion of the transposed segment, analogous to the 
transposed h-set. 


Definition 2.17. Given a segment S between two transversal sections Ejn and Sout 
we define the transposed segment between Sout and Ein by setting: 

15^1 :=| 5 |, 
u(S'^) := s(S), 

^ ^ ' 35 

s(S^) := u(S), 

csT := o o cs; 


where o : x x K —)• x x K, o(p, q, r) 

Observe that 

{s^r = s\ 


^s^)+ = s-. 


{q,p,l-r). 


(36) 


Proposition 2.18. Let S be an isolating segment between transversal sections Tun 
and Yjout for x = f{x). Then S'^ is an isolating segment between Hout and Tim for 
the inverted vector field x = —f{x). The sections Tout and Tm are transversal for 
the inverted vector field. 


2.5. Isolating segments imply coverings. Given a segment S between transver¬ 
sal sections Sin and Eout there is a natural structure of h-sets defined on the faces 
given by intersections Ein H [S’! and Eout O [S'!. 

Definition 2.19. We define the h-sets: 

• C Ein {the front face), 

• .Xs,out c Eout {the rear face), 
as follows: 

• u{Xs, in) = u{Xs, out) ■■= u{S) and s(Xs 4 „) = s(Xs,out) := s{S); 

• |Xs,in| = Ein n |5| and |Xs,out| = Eout n |5|; 

• cxs.in := (7’'„,7r,) ocs|^.^ and := (tt^, tt^) o . 

Definition 2.20. Let S be an isolating segment. We define the exit map Es : 
|-^S,in| —t 5 U |Xs_out| by 

Es{x) = p{te,x), te = uim{t U {0} : :p(t,cc) S 5“ U |Xs,out|} , (37) 


^The sets are manifolds with boundary, so by normals at the boundary points we mean 
normals to any smooth extension of to a manifold without boundary. 
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and the persistent set by 

S° := {x e iXs.inI : Esix) S |Xs,out|}. (38) 

Remark 2.21. From (SI), (S2), (S3) it follows that the function Es is well-defined 
and a homeomorphism onto its image. 


Theorem 2.22. Let S be an isolating segment between transversal sections 
and Yiout- Define V := {x G : 3t > 0 : lp{t,x) G Sou*} and a Poincare map 
P :V ^ Tiout as in equation (29). Then 

• V^$; 

• there exists a dijfeomorphism R : Sm -G such that we have a covering 
relation 


^S,in 


Xs,c 


and 

• it holds that 


P{x) = R{x) Vx G S' 


0 . 


(39) 

(40) 


5° = {x e \Xs,^n\ : R{x) G (41) 

In particular, for every x G such that R{x) G |Xs_o„t| the part of the 

trajectory between x and P(x) = R(x) is contained in [S']. 

The intuition behind this theorem is portrayed in Figure]^ 


Proof. To make the formulas clearer, without any loss of generality we assume that 
cs = idRN. Define 

gix) = ( 42 ) 

where 

/ e^TTu{x) \ 

(fi(t,x)= e~*TTs{x) (43) 

V t + TT^ix) j 

is a global flow on x x M. 

We have 

(V ||7r„(ai)|| ,g(a;)) = ^ \\nufiit,x)\\^^^^ = ||7r„(a;)||, (44) 

(V ||7r,,(ai)|| ,g(a;)) = ^ \\n,(p{t,x)\\^^^^ = - ||7r,,(a;)||, (45) 

for all a; ^ 0 and 

{VTT^{x),g{x)) = = 1, (46) 

for all X G 

Let 17 be a bounded, open neighborhood of [S'], small enough so that the following 
conditions are satisfied: 


(V || 7 r„(a;)|| ,f{x)) > const > 0 Vx S 1/ : || 7 r„(a:)|| > 1, (47) 

(V || 7 rs(a;)|| ,f{x)) < const <0\/x GU : || 7 rs(a:)|| > 1, (48) 

(V 7 r^(a;), /(cc)) > const > QMx GU. (49) 

Let 77 : [0,1] be a C°° function equal to 1 on [S'! and equal to 0 on W^\U. 

Put /(x) = r]{x)f{x) -t (1 — g{x))g{x). Denote by ip{t, x) the local flow generated by 

X = /(x), X G M.^. (50) 
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Figure 3. Schematic drawings of covering relations in a three- 
dimensional isolating segment S with u{S) = s{S) = 1. 


Since = ip{t,x), x G R^\17 and U is bounded, ip is also a global flow. 

From (461 and (49) we have: 


(V7r^(x),f(x)) = r/(x)(V7r^(x),f(x)) + (1 - ?7(x))(V7r^(x), g(x)) 
> const > 0 

for all X G ■ Therefore, the Poincare map 

Pf : X X {0} ^ x x {!}, 

is a well-defined diffeomorphism. We set R := Pj. 


(51) 


(52) 


First, we will prove that i?-covers Xg^out, cf. Figures 3a and 3b In what 
is below we identify the spaces x x {0} and x {1} with 
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two copies x wherever necessary, by projecting/embedding the first 

u{S) + s{S) coordinates. 

We need a homotopy of i? to a linear map. Consider the parameterized family of 
vector fields 

/«(2;) = (1- 0 /( 2 ;)+^g(a;), (53) 

where ^ G [0,1]. By the same reasoning as with / each of these vector fields generates 
a global flow and induces an associated Poincare map 

Pf^ : X X {0} ^ x x {!}. (54) 


We define a continuous homotopy of maps h : [0,1] x x —)• x 

h{i,-)--=Ph„ [ 0 , 1 / 2 ], 


MC,-):= 


;idm 


*(S) 


(2 - 2^)e-i ida 


,(S) 


Ce [1/2,1]. 


(55) 


Indeed, the homotopy agrees at 1/2. Moreover, h{0, ■) = R and h{l, •) is a linear map 
satisfying the requirements given by (14). Since it is also clear that (12) and (13) 
hold for ^ € [1/2,1], we proceed to check these two conditions on the other half of 
the interval. 

Denote by the family of global flows generated by cfc = f^{x). From (47) 
and (44), for ^ S [0,1] and x : [[^^(a;)]] > 1 we get 


= (Vll7r„(a;)ll,/^(a;)) 

= (1 -0(V ll7r„(a;)ll ,/(a;)) + ^(V ll7r„(a;)ll ,5(2;)) 

= (1 -0^(a;)(V ll7r,,(ai)ll ,/(x)) + (1 - r]{x) + ^r]{x)){V \\Truix)\\, g{x)) 
> const > 0. 


Therefore, [[^^(a;)]] = 1 implies [ [7r„(P^^(a:))] [ > 1 for all ^ S [0,1] and proves 

< const <0, X : lj7rs(a:)ll > 1, (57) 


(12). By a mirror argument, from ( |48[ ) and (45) we obtain 
d 


dt 


{t,x))\ 


|t=0 


hence ||7rs(P/5(a;))] ] = 1 implies [[^^(a;)]] > 1 for all ^ S [0,1]. This proves (13). 

We are left to prove (40) and (41). Let us start with the latter. 

Observe, that which proves the “C” inclusion. For the other one 

we proceed as follows. Since (p^ = tp, from (56) and (57) we obtain the forward 
invariance of the sets {x G : ll7r„a;ll > 1}, {a; G jj7rs(x)ll < 1} under ip. 

Therefore, for x G l^s,in| such that R{x) G [Xs^out] we have 


(58) 


\\TTup{t,x)\ \ <1, Vt > 0 : -K^p{t,x) < 1, 

ll7rs(/3(t,a;)ll <1, Vt > 0 : > 0. 

As a consequence the part of the trajectory between x and Pj{x) is wholly 

contained in [S'], where the vector field / is equal to /, hence x G 5°. By the same 
argument (|4l]) implies (40). □ 


One can also prove an analogous backcovering lemma, which is superfluous in 
the context of our applications and therefore will be omitted. 

We are now ready to state a prototypical theorem on how a sequence of covering 
relations and isolating segments forces the existence of periodic orbits for a given 
differential equation. We do not use it though for the FitzHugh-Nagumo equations; 
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a more refined version, which includes a switch between the slow and the fast 
dynamics - Theorem |2.30| - is applied instead. 

Theorem 2.23. Let x = f{x), x G be given by a smooth vector field. Assume 
that there exists a sequence of transversal sections {Sijfhg, /c € N and a sequence 
of h-sets 

X = {X,: |X,| C i = 0,...,k}, (59) 

such that for each two consecutive h-sets Xi-i, Xi G X we have one of the following: 
• there exists a Poincare map Pi : > Si with f2i_i C Si_i and an 

integer Wi G Z* such that 


X^-i <= 


X,,. 


(60) 


• there exists an isolating segment Si between Ei_i and Si such that Xg^ i^ = 
Xi—\ and Xg.^^out — Xi. 

Then, there exists a solution x{t) of the differential equation passing consecutively 
through the interiors of all Xi’s. Moreover: 

• whenever Xi-i and Xi are connected by an isolating segment, the solution 
passes through 5°; 

• if Xq = Xi- the solution x{t) can be chosen to be periodic. 


Proof. By applying Theorem 2.22| we get a chain of covering relations 


X, 


0 ^ 




32,^2 


x. 


9k,Wk 


Xk, 


(61) 


where gi = Pi or gi = Ri, Ri being the diffeomorphism given by Theorem |2.22| 
associated with the segment Si (then Wi = ±1). From Theorem 2.8 there exists a 
sequence {xi : Xi G int |Xi|,z = 1,... ,k} such that gfixi-i) = Xi and we can choose 
a^o = Xk whenever Xq = Xk. 

Suppose that for certain z’s we have gi = Ri. Since Xi-i G |-^i-i| and Rfixi-i) = 
implies that Xi-\ G S'°_i and Rfixi-i) = Pfixi-i), Pi : 


Xi G 


2.22 


\Xi\, Theorem 

—> Si being a Poincare map defined on a subset of Si_i. This proves that this 


orbit is an orbit of a full sequence of Poincare maps, hence a real trajectory for the 
flow. Furthermore, it is a periodic trajectory if xq = Xk (notice that it cannot be an 
equilibrium as the vector field on transversal sections cannot equal 0). □ 


Corollary 2.24. For an isolating segment S the set S^ is nonempty. 


2.5.1. Additional coverings within an isolating segment - the “fast-slow switch”. 
Let us first explain the ideas behind this subsection without formality. Consider a 
three-dimensional isolating segment S with one exit and one entry direction let us 
write ^5 lu, Xg ru for the two connected components of the exit set S~ a “left exit” 
and a “right exit” one, respectively. Each of them lies within a level set given by 
fixing the exit direction level to =Fl. They can be equipped with an h-set structure 
with one exit and one entry direction by setting the entry direction of the segment 
as the entry one and the central direction of the segment as the exit one. 

If we now consider the function Es, which maps each point of the front face 
to the point of dS where the trajectory leaves S, then its image will give a 
similar alignment as in Lemma 2.9 see Figures 3a and 3c The left/right exit edges 
of Xggn remain stationary and coincide with the left exit edges of Xg lu, ^s,ru, so 
to get an actual covering one needs to constrict the h-sets in the image in the exit 
direction by a small factor. 

For the two connected components of S'+ - the “left/right entrance” h-sets 
Xs^rs one needs to fix the entry direction height in the segment coordinates so the 
central direction of the segment takes its role, while the exit direction of the segment 
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induces the exit direction for the h-set. Then one can prove similar theorems with 
backcovering relations, by reversing the vector field. 

In the context of the FitzHugh-Nagumo model such relations allow us to describe 
the passage between slow and fast dynamics where the periodic orbit detaches 
from the slow manifold and starts following a heteroclinic connection of the fast 
subsystem. With an eye on this application we will state the subsequent results 


for a range of dimension combinations which allows an easy proof by Lemma 2.9 


We suspect similar theorems hold for all dimension combinations, and it will be a 
subject of further studies to formulate adequate proofs. 

Definition 2.25. Let S' be a segment with u{S) = 1 and s{S) = s. We define the 
h-sets: 

• -^S.iu C Cg^({— 1 } X X K) {the left exit face), 

• -^S.ru C Cg^({l} X R® X M) {the right exit face) 
as follows: 

• u{Xs,lu) = u{Xsjn) ■= 1 and s{Xs,in) = s(Xs_ru) := s; 

• we set 

l^s.iul := Csi({-l}xK®xR)n|S|, 

|^S,ru|:= Csi({l}xK®xK)n|S|; 

• we identify {t 1 } x M® x M with IR®+^ and then set 


(62) 


CAs.iu := Pu o cs| 


:= Pu 

where pu{p,q,r) = {2r -l,q). 


({ —1} XM) 
XR) ’ 


(63) 


In the above definition the role of pu is to change the order of coordinates, as the 
third center variable in S becomes an exit variable in Xsju and ^s,ru- 

Definition 2.26. Let S be a segment with u{S) = u and s{S) = 1. We define the 
h-sets: 

• C Cg^(K“ X {—1} X R) {the left entrance face), 

• C Cg^(R“ X {1} X R) {the right entrance face) 
as follows: 

• u{Xs,\u) = u{Xs,rn) ■= U and s(Xsjs) = u{Xs,rs) ■= 1; 

• we set 


we identify . 


|Xs,i,|:=Cg^(R“x{-l}xR)n|S|, 

|Xs,,,|:=Csi(R“x{l}xR)n|S|; 

X {=f 1} X K with R"+\ then set 


CAs.i. := Ps ocs|^-i( 


CXs 




‘x{-l}xIE 

‘x{l}xR) 


(64) 


(65) 


where ps{p, q, r) = {p, 2 r - 1 ). 


The role of ps is to change the center variable in S to an entry variable in the 
h-sets Xgjs and 

Definition 2.27. Let X be an h-set with u{X) = u and s(X) = s and let 6 > 0. 
We define: 

• the (5-constricted in the exit direction h-set 

• the (5-constricted in the entry direction h-set 
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by setting: 


where Vuc,Vs. 


— Vuc ^ 

= 'Use ^ Cx 5 

-A K“ X K® and: 

Vnc{p,q) = ((1 + 5)p,q), 
Vsc{p,q) = {p,(l + 6)q). 


( 66 ) 


(67) 


Geometrically, (5-constriction shortens the h-set by a factor 1/(1 -|- (5) in the 
exit/entry direction. Our notation uc, sc stands for constricted in the “unsta¬ 
ble” / “stable” (i.e. exit/entry) direction. 


Theorem 2.28. Let S be an isolating segment between transversal sections 
with u{S) = 1 and s{S) = s. We have the following covering relations: 


Xs,i'. 

Xs,K 

Es. 


Es , 


Es ^ 


-^d,uc 


X 


6,uc 

S,ru'> 


and 


( 68 ) 


for all S > 0. 

Proof. We will only prove Xs^m other case is analogous. The idea 

of the proof should become immediately clear by looking at Figure]^ We embed 
the codomain of Es in a folded a folded hyperplane consisting of three parts: 

• The “upper part” := ({1} x K® x (—oo, 1]); 

• the “middle part” P^'su ■— ([“^i 1] x M® x {1}); 

• the “lower part” := ({—1} x M® x (—oo, 1]). 

It can be regarded as a piecewise smooth section homeomorphic to IR®+^, transversal 
in the sense that there exist smooth extensions of its smooth pieces Sg„, E™^, Eg„ 
to manifolds without boundary which are transversal sections for the vector held. 

We equip Tjs,u with a coordinate system which is given by any homeomorphic 
extension of coordinates given on E^^ by cxg to all 'Es,u - we denote this 

extension by ^. 


The plan is to use Lemma 2.9 and prove conditions that give the same topological 
alignment as needed for a covering relation. 

Recall, that by Es^c we denote the exit map expressed in local coordinates of the 
h-set Xgjn and the section Es_„. In the Es^„ coordinates the support of is 
a product of two balls yqij] x Bs{0, 1). To be in formal agreement with the 
dehnition of the support we would need to stretch out the hrst ball to [—1,1] but it 
is clear that assumptions of Lemma |2.9| are given by geometrical conditions which 
persist under such rescaling. Therefore we omit this transformation to keep the 
notation simple. 

By definition of Eg^^ we have 

l^s.inl n Es,, = (69) 

Coupled with the coordinate system we have chosen on 


hence Es\x- = id 
Ss,« we get 




Es,e ({-1} X R,(0,1)) = {-1} X 5,(0,1), 
EuEs,c ({1} X 5,(0,1)) > 1. 


( 70 ) 
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This, after the aforementioned rescaling of the exit coordinate, implies Condition 


(C2) in Lemma 2.9 


Condition (Cl) follows easily. From (691 and (S3a) we have 

7T. {EsAXs,in,c) n [-1/(1 + 6), 1/(1 + 5)] X W) C B,(0, 1), 
since we need non-zero positive time to reach 'Eis,u- 


(71) 

□ 


If we consider the exit map Egr for the reversed flow x = —f{x) in a transposed 
segment we obtain the following theorem. 


Theorem 2.29. Let S be an isolating segment between transversal sections and 
Eout with u{S) = u and s{S) = 1. We have the following covering relations 


. sr -yS,SC 

^S,out 

^ , sr -yS,SC 

^S,out < ^S,rs 

for all S > 0. 


(72) 


Theorem 2.30. Let x = f{x), a; € be given by a smooth vector field. Assume 
that there exists a sequence of transversal sections {Si}lCQ, k gN, and sequence of 
h-sets 

T = {X, : uiX,) = six,) = 1, i = 0,..., fc} (73) 

such that for each two consecutive h-sets Xi-i, Xi G X we have one of the following: 

• Xi-i C Si- 1 , Xi C Si and there exists a Poincare map Pi : Oi_i -G Si 
with rti-i C Si_i and an integer wt G Z* such that 

X,_i X„ (74) 

• there exists an isolating segment Si between Si_i and Si such that Xs,,i„ = 
Xi—i and Xg^,Qy,i — Xi, 

• there exists an isolating segment Si between Si_i and Si such that Xg.,in = 
Xi_i and either Xg.jy, = Xi or Xg,,ru = Xi; 

• there exists an isolating segment Si between Si_i and Si such that Xg,,out = 
Xi and either Xg,,is = Xi-i or Xg,,rs = -^i-i- 

Then there exists a solution x(t) of the differential equation passing consecutively 
through the interiors of all Xi’s. Moreover: 

• whenever Xi-i and Xi are connected by an isolating segment as its front 
and rear faces, respectively, the solution passes through S^ 

• if Xq = Xk the solution x(t) can be chosen to be periodic. 


Proof. First, we replace all the h-sets Xi of the form Xg.^iu, by the constricted 

versions -^ 5 /))/, -^g/ru ^he h-sets of the form X 5 .JS, ^Si,rs by -^ 5 / 1 / Xgl^s- 
Let us denote the new h-sets by Xi. The replacement procedure is done one by one. 
Each time an h-set Xi needs to be replaced we choose (5i > 0 small enough, such 
that 

( 1 .) any covering relation Xi was involved in is preserved for Xi, 

( 2 .) any isolating segment that was built including Xi as either the front or the 
rear face can be reconstructed as an isolating segment Si/Si+i with the face 

It is intuitively clear that both should hold for a sufficiently small perturbation. To 
show ( 1 .) it is enough to observe that a covering relation is a (7°-open condition 
with respect to homeomorphisms defining the h-sets and persists after constricting 
one (or both) h-sets with <5 small enough. The proof of such proposition would be 
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almost the same as the proof of Theorem 13 in stating stability of covering 
relations under (7° perturbations, and therefore we omit it. 

For (2.) the segment Si is constructed so that cs^ is 0(i5i)-close in the norm 
to Cg.. We omit the details; describing the construction by precise formulas would 
introduce a lot of unnecessary notation. It is easy to see that for Si small enough 
the conditions (S1)-(S3) (or their counterparts) will still hold. 

We apply Theorems |2.22[ |2.28[ |2.29| to get and a chain of covering relations 


X, 


0 ^ 


Xi 






9k:'^k. -G 
< ■ • ■ > Afc, 


(75) 


where for each gi we have one of the following: 

• gi = Pi, 

• gi = Ri, Ri given by Theorem |2.22[ 

• 9i = ESi 


^ 9i — ^s '^' 

From here, the proof continues in the same way as the proof of Theorem |2.23[ We 
obtain a sequence of points : Xi G intXi,i = 1 ,..., fc} such that 9 i(xi-i) = Xi 
and we can choose xg = Xk whenever Xg = Xk- By the same argument as in 
Theorem |2. 23 the sequence lies on a true trajectory of the flow; the trajectory is 
periodic if xg = Xk- □ 


We note that the formulation of Theorem |2.30 is not aimed at full generality. 
By using only Theorem |2.28| or |2.29| one can produce similar theorems when one 
direction is expanding and arbitrary number of directions are contracting or vice 


versa. 


3. A MODEL EXAMPLE 


The purpose of this section is to discuss a model example for the construction of 
a closed chain of covering relations and isolating segments in the FitzHugh-Nagumo 
equations. Contents of this section are by no means necessary to prove the main 
Theorems onioi Instead, their purpose is to predict that the computer 
assisted proof based on Theorem |2.30| will succeed, so the sequence of coverings is 
found not just by pure luck, but is backed up by the structure of the singularly 
perturbed system. 

Our example will be a fast-slow system of the form: 


y = ^9{x,y,P), 


(76) 


where a: G y G M, /, g are smooth functions of {x, y, e) and 0 < e <C 1 is the small 
parameter. We will also write x = {xi,X 2 ) to denote the respective fast coordinates. 
By the slow manifold of (76) we will mean the set {{x,y) : f{x,y,0) = 0}, and by 
the fast subsystem a two-dimensional system 


i = f{x,y,0), 


(77) 


given by fixing ?/ as a parameter. 

We make the following assumptions: 

(PI) we have two branches of the slow manifold A±i, that coincide with {0} x 
{1} X K and {0} x { — 1} x K, respectiveljQ Both are hyperbolic with 


^One expects that the branches would actually connect with each other, but we bear in mind 
that this is a model example and the fold points are of no interest to us. 
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one expanding and one contracting direction, and the vector field in their 
neighborhoods U±i is of the following form: 

f{x,y,0) = A±i{y){xT [0,1]’^) + h±i(x,y), (78) 

3eo > 0 : 0 < i5±i < ±g{x, y, e) < , {x, y) S t/±i, e G (0, cq]- (79) 


The functions A±i, h±i are assumed to be smooth and to have the following 
properties 


A±i{y) = 


A«,±i(y) 0 

0 K,±i{y) 


—<^±i < As_±i(i/) < —(5±i < 0 < 5±i < A„_±i(y) < 


h±i{x,y) x^{Q,±i) 

Ik-(0,±1)11- 


(80) 

(81) 

(82) 


The values (5±i > 0 are some constant bounds, which in particular do not 
depend on neither e nor y. 

(P2) For the parameterized family of the fast subsystems we have two parameters 
y*,y*, without loss of generality assumed to be equal to =Fl, for which 
there exists a transversal heteroclinic connection between the equilibria 
(0, —1) and (0,1) in the first case, and (0,1) and (0, —1) in the second. That 
means: given any two one-dimensional transversal sections for the 

fast subsystems for y = ±1 which have a nonempty, transversal intersection 
with the heteroclinic orbits, the maps given by 


'I'ii: y —t w^±LE/,±i(y) — ^Ti±/,±ik) 


(83) 


have zeroes and a non-zero derivative at y = ±1. 

Here ^{y) and ^{y) denote the first intersections between the 
appropriate branche^ of the unstable/stable manifolds of the equilibria 
(0,±1) with a given section S in the section coordinates. 

(P3) Denote the points (0, —1, —1), (0,1, —1), (0,1,1), (0, —1,1) by Pq,, a gI = 
{dl, respectively and set e = 0. For each a G I there exists 

a neighborhood Va of Fq,, such that if A±i n Va is the part of the slow 
manifold contained in Va, then the part of its unstable manifold contained 
in Va coincides with the plane M x {±1} x M, and the part of the stable 
manifold contained in Va - with the plane {0} x K x K. Without loss of 
generality we can have UaGi^“ (^-i U C/i). 


Assumptions that provide us with straightened coordinates are used mostly 
to simplify the exposition. It is our impression that the Fenichel theory, and in 
particular the Fenichel normal form around the slow manifold are well-suited for 
verifying such conditions, see 


12 


Theorem 3.1. Under assumptions (P1)-(P3), there exists an eq > 0 and six sets 
forming isolating segments for ( [7^ for e G (0,eo].' 

• Su, Sd - two “long” isolating segments positioned around the branches A±i 
of the slow manifold; 

• Sa, OL G T - four short “corner” isolating segments, each containing the 
respective point F„; 


^To not complicate further the notation, we make an implicit assumption that only one pair of 
branches cross in each of the two subsystems and only refer to them. 

^ Index letters in X stand for up/down and left/right and refer to positions of the points in the 
(y,X 2 ) plane, see Figure]^ 
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Figure 4. Isolating segments and Poincare maps in the model 
example for the periodic orbit, the sequence of h-sets plotted in 
red. 


along with the associated transversal sections of the form 

u{Sdl) = s{Sdl) = u{Sdr) = s{Sdr) 

— — s(^Sul^ = U^Sur') — ^{Sur') 

= u{Su) = s{Su) = u{Sd) = s{Sd) 


with 


(84) 


= 1 . 


Moreover, for the h-sets defined by isolating segments we have 


^Su.,OUt — ^Sur,i‘fT--; 

(85) 

^SdriOUt - ^Sd,'i'TI'‘! 

(86) 

^Sd,out — ^Sdi.’i-ri: 

(87) 

^Sui,out — ^Su,inj 

(88) 


and the collection 

Xfhn,p = 

^outl (89) 

XSd,inj XSd ,OUt-) ^Sdi,ru, ^S„i,lsi 


satisfies assumptions of Theorem 2.3(J[ for e G (0, cq] . In particular we have the 
following covering relations among the h-sets not connected by an isolating segment: 


USdl.ru 

Pl , -y- 

(90) 

^S„.r,lu ■ 

Pr. 

-^ -^Sdr.rs, 

(91) 


where P* are Poincare maps between the respective h-sets and transversal sections 
containing the next h-set. 

As a consequence there exists a periodic solution of the system for these parameter 
values. 


The conclusion of the theorem is portrayed in Figure]^ We break the proof into 
two parts, first we prove the existence of the corner isolating segments and coverings 
as a separate lemma. 
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Lemma 3.2. Consider the system (76). For e S (0, eo], eo > 0 small there exist 
two transversal sections of the form 


Sr := {{xi,X 2 ,y) : X 2 = 1 - El} C V^i C Kz, 
Sr := {{xi,X2,y) : X 2 = -!+£«} H Vdr C Vdr, 


(92) 


Vui,Vdr being neighborhoods ofTui and and four isolating segments Sdi, Sui, 
Sur, Sdr as specified in Theorem\3.1\ sueh that 


C Sr, 
C Sr, 


(93) 


and there are coverings 


^Sdl,ru 

Pl , Y' 

(94) 

^S„r,lu '■ 

Pr. V 

-^ XSd,,rs- 

(95) 


Moreover, the seetions S* and the segments 5"* are e-independent, and given a 
maximal diameter diarrimo^: > 0 they can be chosen so that 


diam(S'*) < diammax- 


(96) 


Proof. We focus first on the “left” part of the picture, since all arguments for the 
“right” part are symmetric and independent. Without loss of generality we can 
assume the crossing of the unstable and stable manifolds near the point F^; occurs 
for X 2 — 1 negative and take Sl > 0. For £r and e small enough, condition (PI) 
implies that the linear part of the vector held dominates the higher order terms 
h±i, so after having set a sufhciently small neighborhood Vui the section Er is 
transversal. 

The construction of the isolating segments Sui,Sdi is also enabled by (PI). 
Because we already work in straightened coordinates, their supports can be chosen 
to be of the form: 

\Sdi\ = [—El,El] X [—1 — Ss^di, —1 + x [—1 — 6u,di,—i + Su,di], 

I ^ul I [ , Sy^^. 1 , 1 ] X [1 £r , 1 -)- £r] X [ 1 dgW, 1 “f , (1^^) 

u{Sdi) = s{Sdi) = u{Sui) = s{Sdi) = 1, 


where the constants Sg^di, bu,di,bs,ui, bu,ui will be hxed later in the proof. The changes 
of coordinates , cs„, are dehned as a translation of the cuboids to the origin of 
the coordinate system composed with rescaling to [—1,1]^ x [0,1]. We label the 
hrst coordinate as exit, second as entry, third central. Again, if el, Sg^di, du,di, du,ui, 
Ss^ui are small, then the linear part of the vector field dominates the nonlinear part 
and conditions (S2b), (S3b) are satisfied for e = 0 and for e > 0 small. Since our 
change of coordinates is of the form as in (33), for e > 0 small (Sla) follows from 
the inequalities 


We can now move on to proving the covering relation (94). The supports of the 
h-sets are of the form: 


I^S'di,ru| X [ 1 bg dl, f T bg dl] X [ 1 di, dl, 1 T ^u,cZz]; 

I ,1s I — [ du,ul, du.ul] X {1 Eif\ X [ 1 1 T f^Sj'uz]- 


(98) 


In X^ji^ru the X 2 variable takes the role of the entry variable and y takes the role of 
the exit one; in the variable xi is exit and y is entry. 

Since covering relations are robust with respect to perturbations of the vector 
field (see Theorem 13 in 39 ) it is enough to show them for e = 0. From (P2) and 
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(P3) we know that 




A 

dy 






( 99 ) 


and without loss of generality let us assume that is increasing in the neighbor¬ 

hood of the point (e^, —1, —1). That allows us to define two h-sets with a covering 
relation between them. The procedure is as follows: 




fix some Ss^ui > 0. 

To comply with the covering condition (Cl) from Lemma 
and Su.di > 0 so that 


2.9 


choose Ss^di > 0 


nyPL{\Xs^,,is\) c (-1 - 6s,ui,-l + Ss,ui). (100) 


Now, provided 6s^di and 5u,di were chosen small enough, 
exists Sul > 0 such that 


from (99) there 


'^x\Pl ({s:^} X [ 1 ds^dii 1 T ^ { f du^di}') ^ £ui ^ 0, 

'^x^Pl ({s:^} X [ 1 Ss^di: 1 T ^ { 1 T ^ £ui ^ 0. 


( 101 ) 


• To fulfill (C2) it is enough to choose Su,ui < £ui- 

It is clear that we can choose Sl small enough and then perform the procedure 
above with 6’s small in a way, that the diameter constriction (961 is satisfied. 

The same procedure is repeated for the isolating segments Sdr] we will only 
introduce the notation for these segments, as they will be used later in the main 
part of the proof of Theorem |3.1| Similarly to the left side segments we define them 
by giving the cuboid supports 


I * 5 ''ur I — [ £Ri£r\ ^ [f 1 ^ [f 1 

I — [ du,dri du^dr] ^ [ 1 £R^ ITS/j] X [1 (102) 


and the coordinate changes are again simple translations and rescalings 

to [—1,1]^ X [0,1], so the first variable in the supports is the exit one and the second 
is entry. 

The supports of h-sets of interest are as follows: 

|^Su,.,lu| — { £r} ^ [1 1 T ^ [1 du^ur-, 1 T du^ur\^ (103) 

|^5d,.,rs| — [ du^drj ^ { 1 T £R^ ^ [ 1 ^s,dr j 1 “t” 

We will not go into details of determining Su,ur, dg^ur, du,dr, ds,dr and er - the 
procedure is exactly the same as for the left side segments. The variable y is the 
exit variable and X 2 is the entry variable in as for xi is the exit 

one and y is entry. 

By taking the minimum of all upper bounds on e’s throughout this lemma we 
obtain eg and the proof is complete. □ 


We can now return to proving Theorem |3.1| We import all the notation from 
the proof of the Lemma |3.2| and in particular assume that the isolating segments 
Sdi^ Suij Sun Sdr and the respective h-sets can be chosen to be of the form given 
therein. 


Proof of Theorem 3.1 
diameter diammax 


From Lemma 3.2 for any given maximal corner segment 


> 0 we obtain a bound eg on e’s and four isolating segments 
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Sdi, Sul, Sur, Sdr containing the respective points Fq with covering relations between 
their respective faces. We set diammax small enough to have 

f{x,y,0) Ki Ai{y){x - [0,1]'^), (t, y) e conv(|S'ui| U (104) 

f{x,y,0) A_i{y){x+[0,lf'), (a:,y) S conv(|S'di| U (105) 

so the higher order terms h can be assumed negligible when checking isolation 
inequalities in these neighborhoods. 

Given our four corner isolating segments we are left with construction of two 
isolating segments Su and Sd which connect the pairs Sui, Sur and Sdr, Sdi, respec¬ 
tively. We will only construct Su, the case of Sd is analogous. The strategy is to 
first connect the pairs by segments, then, if necessary, decrease cq to some smaller 
eo to obtain isolation. 

We introduce the following notation for rectangular sets around the upper branch 
of the slow manifold: 


LuiSu, Ss,y) ■■= [-5„, 5u] X [1 - 5s, 1 + X {y}. 


(106) 


We set 


Otj —1 -|- 5s 

5u ■— 1 5u.u 


,ul 1 


(107) 


and we can assume that Uu < bu- Now, we can define Su as a cuboid stretching 
from to as follows. For the support we put 

\Su\ := 

[J ((1 — -I-(1 — (1 — + ■f^u) • 

«e[ 0 T] 

We also set u{Su) = s{Su) ■= 1- There is no need for description of cs^ by precise 
formulas, so we only mention that it is a composition of 

• a diffeomorphism which rescales each fiber (•,!/), given by fixing y S 
[au,bu], to [-1,1] X [-1,1], 

• a rescaling in the central, y direction from [a„, 6„] to [uu, -b 1], 

• a translation to the origin of the coordinate system. 

As with the corner segments, xi is labeled as the exit direction, X 2 as entry, and 
y as the central direction. Then one sees that equalities (85) and (88) are true. 
Condition (Sla) is a consequence of inequalities (79) for small e , as the change of 
variables cs^ in the central direction takes the form (33). The upper bound for e’s 
given by eo may need to be decreased at this step. 

It remains to check (S2b) and (S3b) and for that purpose we may need to further 
reduce eg. Normals to S'” pointing outward of |S„| are given by 

5u.ur 5u ul 


n-{x,y) = sgna;i, 0 , - 


bu Rk 


From (104), ( 109[ ) and (80), (81) for {x,y) G Su we have 
{Uxg),n-){x,y,e) « A„,i(y)|a:i| - eg{x,y,e) 


5u,.ur 5, 


u,ul 


^7/, ' 


> 5i|a::i| - 




(109) 


( 110 ) 


5i bu du 

and the right-hand side is greater than 0 for e G (0, eg], eg small enough, see Figurej^ 
This proves (S2b). 

Verifying (S3b) goes along the same lines, the expression for outward normals is 

5s.ur 5s ul 


n+ix,y) = 0 ,sgn(a :2 - 1 ),-- 


ba, a,. 


( 111 ) 
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y 


dt 


Figure 5. Isolation in segments around the slow manifold for 
small e. The fast component of the vector field dominates the slow 
one and offsets the influence of the slope on isolation inequalities. 

and one readily checks that 



( 112 ) 


decreasing eo if necessary. 

The only difference in the construction of Sd is that the recipe for has to include 
a flip in the y direction so we can have ^Sd,in = ^Sd^,out and 
By taking minimum of all upper bounds for eg throughout the proof we obtain the 


desired eq. 


□ 


4. Computer assisted proofs of Theorems oo AND |1.3| 

Most of the numerical values in this section are given as approximations with 8 
significant digits. An exception to that are the equation parameters, which are exact. 
Therefore, computations that are described below are not actually rigorous, but 
the programs execute rigorous computations for values close to the ones provided. 
Actual values in the program used for rigorous computations are intervals with 
double precision endpoints - we decided that writing their binary representations 
would obscure the exposition. If needed, exact values can always be retrieved by the 
reader from the programs. If interval is very narrow and used to represent only one 
particular value, such as a coordinate of a point, we just write a single value instead. 

Rigorous and nonrigorous integration, computation of enclosures of Poincare 
maps defined between affine sections and their derivatives, linear interval algebra 
and interval arithmetics is handled by routines from the CAPD libraries and 
we do not discuss it here. For rigorous integration we used the Taylor integrator 
provided in CAPD. 

The source code executing the proofs is available at the author’s webpage [^. 
Our exposition loosely follows what is performed by our programs. The best way to 
examine the proof in detail is to look into the source code files. For most objects we 
use the same notation in the description as in the source code, however occasionally 
these two differ. In such cases identifying the appropriate variables should be easy 
from the context and from the comments left in the source code files. 

For a given vector object x, by x[i-l] we denote the i-th coordinate of x. We 
will denote the right-hand side of (U) by F. 

We recall that the fast subsystem of (U is given by: 


u = V, 


v' = 0.2{0v — u{u — O.I)(I — u) -I- w). 


( 113 ) 
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Unless otherwise specified, the half-open parameter intervals e S (0, eo] is treated 
in computations by enclosing it in a closed interval [0, cq]. The assumption e ^ 
0 is utilized only in verification of condition (Sla) for isolating segments (see 


Subsection 4.1.2). 


4.1. General remarks. 


4.1.1. H-sets and covering relations. Each h-set X appearing in our program is 
two-dimensional with u{X) = s(X) := 1 and can be identified with a parallelogram 
lying within some affine section. Verification of covering relations is done exclusively 


by means of Lemmas 2.9 2.11 (see also Remark 2.10). The procedure is relatively 


straightforward and has been described in detail in several papers, see for exam¬ 
ple |35| , therefore we do not repeat it here. We only mention that, if needed, the 
procedure may include subdivision of h-sets. This reduces the wrapping effect, but 
greatly increases runtimes (note that wrapping is already significantly reduced by 
use of the Lohner algorithm within the CAPD integration routines). Given an h-set 
X we want to integrate with subdivision, we introduce an integer parameter div. 
It indicates into how many equal intervals we divide the set in each direction. For 
example, setting div = 20 means that we integrate 20 pieces of and 

400 pieces of |V| to evaluate the image of the Poincare map. In the outlines of 
our proofs we will indicate the values of div to emphasize which parts of the proof 
involved time-consuming computations. 

4.1.2. Segments. Our segments are cuboids placed along the slow manifold Cg so 
that a part of it belonging to the singular orbit is enclosed by them. For each 


segment S we have u{S) = s{S) = 1. All of the segments have the property (33) 


with the slow variable w serving as the central variable. Therefore, to establish 


(Sla) it is enough to show (34), which is equivalent to verifying either u > w or 


u < w for all points of the segment. This in particular allows us to handle half-open 
ranges e G (0, eg] computationally, as at this point we effectively factor out the small 
parameter. 

Confirming (S2b) and (S3b) is simple, as all of the faces lie in affine subspaces. As 
the exit/entry directions we take the approximate directions of the unstable/stable 
bundles of Cg. Similarly as for verification of covering relations, we subdivide the 
sets 5”“, S'^ before evaluating isolation inequalities. The normals are constant 
within a face, the actual beneht is in reduction of wrapping in evaluation of the 
right-hand side of the vector field over a face. 

Our segments are rigid and the stable and unstable bundles of Cg actually slightly 
revolve as we travel along the manifold branches. By using a single segment to 
cover a long piece of the branch we could not expect conditions (S2b), (S3b) to 
hold anymore. Therefore we use sequences of short segments, the position of each 
is well-aligned with the unstable/stable bundles of Cg - we call them chains of 
segments. They are simply sequences of short segments placed one after another, so 
a longer piece of the slow manifold can be covered. We require that each segment Si 
from a chain is an isolating segment and that for each two consecutive segments 5^, 
Si-^-i in the chain the transversal section S^^out containing the face Si^out coincides 
with the section containing 5'i+i.in and there is a covering relation by the 

identity map 

AlSi,out (114) 

In other words, the covering relation is realized purely by the change of coordinates 
-1 


CXs 


o c 


For purposes of checking the assumptions of Theorem 
treat the identity map as a special case of a Poincare map, see Subsection 




2.30 

'n 


we 
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A topic we think is worth exploring, is whether chains of segments are a viable 
alternative to numerical integration in computer assisted proofs for differential 
inclusions arising from evolution PDEs; or of stiff systems where one has a good 
guess for the orbit from a nonrigorous stiff integrator. In future we plan to conduct 
numerical simulations to get more insight on that matter. 

4.1.2.1. Representation of segments. Each segment S in our programs can be repre¬ 
sented by 






two points Front, Rear G K^, serving as approximations of points on Cq, 
a 2x2 real matrix P representing the rotation of the segment around the 
slow manifold (this does not need to be a rotation matrix) - it will contain 
approximate eigenvectors of the linearization of the fast subsystem (113) at 
a selected point from Cq n S, 

four positive numbers a, b, c, d > 0 - the pair (a, b) describes how to stretch 
or narrow the exit and the entry widths of the front face of the segment, 
respectively, and the pair (c,d) does the same for the rear face. 


n 


For a pair of points (a, b) and a 2x2 matrix A we define an auxiliary linear map 

a,b,A : ^ by 


Our segment is then defined by 


■ A ' 


aXu 

0 


bxs 


(115) 


Cg = (1 - a;;_,)(Front -F na,b.p(a;«, a:*)) 

-F a;^(Rear -F nc,d,p(a;„, Xs))- 


(116) 


For such segments one can define their front & rear faces and the left/right 
entrance/exit faces ATs in, As^outi As is, Xs^rsj As,iu, As^ru as in Subsection 2.4 
4.1.2.2. Construction of chains of segments. Our recipe for creating a chain of 
segments along a branch of the slow manifold is as follows. We 

assume we are given two disjoint segments So, Su+i positioned along the slow manifold 
Co that we would like to connect by the chain. Without loss of generality we may 
assume that we are on the upper branch of Cq, so |So| is to the left of |Sh+i| in terms 
of the w coordinate. 

For each segment Si we will use its representation 


(Fronti, Reari, Pi, ai, bi, Ci, di) (117) 

given in Paragraph |4.1.2d] Wherever we mention an identity map id between two 
h-sets, we mean the identity map restricted to the common transversal section. 
Our chain will connect the segments So, Su+i in the sense that 

Aso As„i„, (118) 

As,j,out — As,j^^ in. (119) 


We remark that we connect the faces of two segments as this is what we later do in 
the proof of Theorem 0 but with little changes these could as well be any two 
parallelogram h-sets placed on sections crossing Cq. 

Creating a chain is a sequential process akin to rigorous integration with a fixed 
time step; to construct the segment S^ we need to know the representation of the 
segment Si-i, If 1 < i < A we define the segment Si as follows 




we set Fronti := Reari-i. 

The point Reari is constructed by locating an (approximate) equilibrium of 
the fast subsystem (113) with Newton’s method for 


w := Fronti [2] -F l/N, 


(120) 
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and then embedding it into the 3D space by adding the value of w as the 
third coordinate. 

• Columns of the matrix are set as approximate eigenvectors of linearization 


of (113) at Reari. 
For (ai,bi) we put 


ai := Ci-i/factor, 
bi := factor x di-i 


( 121 ) 


where factor is a real number greater than 1. In our case hardcoding 
factor 1.05 gave good results. 

For (ci,di) we put 


i N-i 

Ci “a-N+i H—u~c:o, 
N N 


di :ubN+i 
N 


N-i 

N 


( 122 ) 


-do. 


For the segment Su we proceed by the same rules with the exception that we set 

Rearu Front[j+i, 

Pn '■= Pn+1, 


(123) 


to comply with (119). 

For such Si we check the conditions (Sla), (S2b), (S3b) and the covering relation 
Then, we proceed to the next segment. 

For N large it is easy to satisfy (S2b), (S3b) for each short segment Si, as each 
Pi approximates the directions of the unstable and stable bundle of Cq. Moreover, 
because Reari_i, Rear^ are close, for each i S {1,... ,N} we have 

PioPi_i«id. (124) 

Thus, for the identity map in the h-sets variables we get 

factor 


ids = cxs. in o c 


-1 


0 


0 

1 

factor 


(125) 


and there are good odds that by use of Lemma |2.9| we can succeed in satisfying the 
conditions (118) and ( |114[ ). 

4.2. Proof of Theorem |1.1[ To deduce the existence of a periodic orbit we check 
the assumptions of Theorem |2.30| Our strategy resembles the one given for the model 
example in Subsectionwhich was portrayed in Figure]^ The main modifications 
are due to numerical reasons: 


• we introduce two additional sections on the trajectories of the fast subsystem 
heteroclinics, in some distance from the corner segments, 

• instead of the “long” segments S'„, Sd we place two chains of segments along 
the slow manifold connecting the corner segments - see Paragraph |4.1.2.2[ 

We divide the parameter range e G (0,1.5 x 10“^] into two subranges (0,10“"^] 
and [10“^, 1.5 x 10“"^]. The procedure is virtually the same for both ranges and the 
only reason for subdivision is that the proof would not succeed for the whole range 
e G (0,1.5 X 10“^] in one go, due to an accumulation of overestimates. Following 
steps are executed by the program for both ranges: 

(1) First, we compute four “corner points” 

GammaDL = (-0.10841296,0,0.025044220) « (Ad(w*),w*), 

GammaUL = (0.97034558,0,0.025044220) « (A„(u;*), w*), 

GammaUR = (0.84174629,0,0.098807631) « {Au{w*),w*), 

GanimaDR= (-0.23701225,0,0.098807631) « {Ad{w*),w*). 
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This computation is nonrigorous; in short we perform a shooting with w pro¬ 


cedure for the fast subsystem (113) from first-order approximations of stable 


and unstable manifolds of the equilibria to an intermediate section; this is an 
approach like in 18 . The matrices given by the approximate eigenvectors of 
the linearization of the fast subsystem at points GaimnaDL, GammaUR, GaimnaUL, 
GaimnaDR are 

1 1 
0.34113340 -0.21913340 


PDL = PUR = 


PUL = PDR = 


1 1 
0.46313340 -0.34113340 


(127) 


respectively. 

(2) We initialize four “corner segments” DLSegment, ULSegment, URSegment and 
DRSegment with data from Tableas described in Paragraph |4. 1.2.1 and check 
that they are isolating segments. For checking the isolation formulas (S2b), 
(S3b) we subdivide enclosures of each of the respective faces of the exit and the 
entrance set into 150^ equal pieces. 


Segment 

Front, Rear 

P 

(a,b) = (c,d) 

DLSegment 

GammaDL ± (0, 0, 0.005) 

PDL 

(0.015,0.012) 

ULSegment 

GammaUL =F (0, 0, 0.005) 

PUL 

(0.01,0.015) 

URSegment 

GammaUR =F (0, 0,0.005) 

PUR 

(0.029,0.019) 

DRSegment 

GammaDR ± (0, 0, 0.005) 

PDR 

(0.007,0.03) 


Table 1. Initialization data for the four corner segments. The pair 
(a, b) determines the exit/entry direction widths of the segments 
and the difference |Front [2] —Rear [2] | the central direction width. 


(3) Unlike in the model example - Lemma 3.2, the two transversal sections we 


integrate to are positioned within some distance from the corner segments. We 
move away from the segments because integration too close to slow manifolds 
poses a numerical problem - the vector field slows too much and the routines 
for verifying transversality fail. More precisely, a section leftSection is placed 
on the integration path between the segments DLSegment, ULSegment and a 
section rightSection on the path between URSegment and DRSegment. 

We define the following Poincare maps: 

• pmDL is the Poincare map from ^DLSegment.ru to leftSection, 

• pmUL is the Poincare map from a subset of leftSection to the afhne 
section containing XuLSegment.is, 

• pmUR is the Poincare map from XuRSegment.iu to rightSection, 

• pmDR is the Poincare map from a subset of rightSection to the affine 
section containing XoRSegment.rs- 

Let now us briefly describe what covering relations we verify. 

We integrate the h-set XnLSegment.ru to leftSection and create an h-set 
midLeftSet C leftSection so that it is pmUL-covered by a small margin 
by -^DLSegment,ru, See Lemma [2^ Then, we integrate the h-set XuLSegment.is 


backward in time to leftSection and verify that midLeftSet pmUL-backcovers 


^ULSegment,Is ■ 

The h-set XuRSegment.iu is integrated to rightSection, and, as in the pre¬ 
vious case, we define an h-set midRightSet C rightSection, such that it is 
pmUR-covered by XuRSegment.iu- Then, we integrate the h-set df DRSegment,rs back¬ 
ward in time to rightSection and verify that midRightSet pmDR-backcovers 

^DRSegment,rs ■ 
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Altogether, we have the following covering relations: 


a:, 


DLSegment,ru 


pmDL 


midLeftSet 


pmUL 


X, 


ULSGgment,ls5 


pmUR . . ^ pmDR ^ 

-^URSegmentjlu Hlid.Righ.'tSG't ~ -A.DR,segmGnt,rs • 


(128) 


Parameter div describing partitioning of h-sets for the rigorous integration was 
set to 20. 

(4) To close the loop, we connect the h-sets XuLSegment.out and ATuRSegment.in by a 
chain of segments UpSegment and XDRSegment.out and ATpLSegm ent.in by a chain 
of segments DownSegment as described in Paragraph |4. 1.2. ^ The number of 
isolating segments in each chain N is set to 80. For verification of the isolation 
conditions (S2b), (S3b) in each chain we partition the enclosures of each of the 
faces of their exit and entrance sets into 110^ equal pieces. 

Many choices of parameters in the program were arbitrary; of most importance 
are the exit/entry/central direction widths of the corner segments given in Table 
For very small e ranges (such as e S (0,10“®], e S (0,10“’’]) various reasonable 
guesses would yield successful proofs, due to the eminent fast-slow structure of the 
equations (cf. Section . However, the range of possibilities would diminish as the 
upper bound on e was increased, and finding values for our final e ranges was a long 
trial-and-error process. This can be explained as follows. For large e’s the periodic 
orbit moves away from the singular orbit, around which we position our sequence of 
segments and h-sets. Moreover, the hyperbolicity of the slow manifold, which plays a 
vital role in the creation of the periodic orbit near the singular limit , decreases as e 
increases. Each time a value of a program parameter was adjusted in an attempt to 
succeed with a particular part of the proof, it was possible that another part would 
fail. For example, increasing the central direction widths of the corner segments 
facilitated the verification of covering relations for the Poincare maps; but too much 
of an increase made isolation checks for the corner segments fail; increasing the 
exit direction widths of ULSegment, DRSegment made the exit direction isolation 
checks (S2b) in segments of UpSegment, DownSegment easier to satisfy but 
had a negative effect on the covering relations; etc. It was particularly difficult to 
simultaneously obtain both isolation for the corner segments and covering relations 
in the fast regime. 

By repeating the process of 

• trying to slightly increase the e range, 

• executing the program with given parameters, 

• should the proof fail, changing the parameters in favor of the inequalities 
which were not fulfilled, at the cost of the ones where we still had some 
freedom, 

we obtained a relatively large range of e € (0,1.5 x 10“^], for which the inequalities 
needed in our assumptions hold by a very small margin. In particular, the right 
bound 1.5 x 10“^ was large enough to include it in a continuation-type proof of 
Theorem |1.2[ performed in reasonable time and without using multiple precision. 

It would certainly be helpful to have that procedure automated. As one can see, 
we are effectively dealing with a constraint satisfaction problem (see [34| ) where 
variables, given by the program parameters have to be chosen to satisfy constraints 
given by inequalities coming from covering relations and isolation conditions. In 
addition, verification of whether constraints are satisfied requires execution of the 
program and is fairly expensive computationally. A suitable algorithm for adjustment 
of parameters to satisfy the constraints would allow to extend the range of the small 
parameter e G (0, eo] even further. We remark that obtaining a large value of eo 
in this proof is crucial for achieving this parameter value with further validated 












32 


A. CZECHOWSKI AND P. ZGLICZYNSKI 


continuation algorithms (like the one in Theorem 1.2) This is due to the fact that 
the period of this unstable orbit is roughly proportional to ^ (see Table which 
makes it virtually impossible to track the orbit by numerical integration methods 
for very small e. 


4.3. Proof of Theorem |1.2[ Our strategy is to check the assumptions of Theo¬ 
rem [2^ for a sequence of h-sets placed along a numerical approximation of an actual 
periodic orbit (not the singular orbit). This can succeed for a very small range of 
e, then we need to recompute our approximation, ending up with a continuation 
procedure. 

We start by generating a numerical approximation vector of 212 points from 
the periodic orbit for e = 0.001 obtained from a nonrigorous continuation with 
MATCONT 17 . From there we perform two continuation procedures, down to 
e = 1.5 X 10^^ and up to e = 0.0015. Each step of the continuation consists 
of a routine proveExistenceOf Orbit performed on equation ([^ with an interval 
currentEpsRange of width incrementSize substituted for e. It can be described 
by the following steps. 


(1) Given an approximation vector initialGuess of pm_count points of the periodic 
orbit obtained from the previous continuation step (in the first step this is 
the MATCONT-precomputed approximation), we initialize a Poincare section 
section! for each of the points initialGuessi by setting the origin of the 
section as the given point and its normal vector as the vector as the difference 
between the current and the next point of the approximation. Then, we refine 
the approximation by a nonrigorous computation of Poincare maps and their 
derivatives and application of Newton’s method to the system of the form (§. 
Note that we set the normal vector to be the difference between the current and 
the next point on the orbit rather than the direction of the vector field, as the 
latter can be misleading close to the strongly hyperbolic slow manifold. Let us 
denote by correctedGuess the Newton-corrected approximation. 

(2) Each section! is equipped with a coordinate system used for the purposes of 
covering by h-sets as described in Subsection |4.1.1| The first column corre¬ 
sponding to the exit direction is obtained by a nonrigorous integration of 
any non-zero normalized vector by the variational equation of ([^ along the 
approximated orbit until it stabilizes; and then propagating it for each i by 
one additional integration loop. Similarly, the second column (corresponding 
to the entry direction) is computed by backward integration of any non-zero 
normalized vector until it stabilizes and further propagation by inverse Poincare 
maps for each i. Then, we project these columns to the orthogonal complement 
of normal!. 

(3) Let pmj^ be a Poincare map from a subset of the section section! to the 
section section!+i mod pm.count- We initialize a sequence of h-sets Xi on sections 
section! by specifying Xq and generating such Xi,..., Apn count-i, so the 

covering relations Xi - -> i € {0,... ,pm_count — 2} hold by a small 

margin. The periodic orbit is strongly hyperbolic and the h-sets will quickly 
grow in the exit direction. Therefore we put an additional upper bound on the 
growth of the exit direction to prevent overestimates coming from integrating 
too large h-sets. For rigorous integration of h-sets the parameter div was set to 
5. 

(4) We check that the following covering relation holds 


A 


P*^ni.count — 1 

~> ^ 0 - 


pm.count —1 


(129) 










PERIODIC SOLUTIONS OF THE FHN EQUATIONS 


33 


This implies the existence of the periodic orbit of for e € currentEpsEange, by 
Theorem 12.81 

(5) We produce a new initialGuess for the next step of continuation by removing 
the points from the approximate orbit where the integration time between 
respective sections is too short and adding them where it is too long. This way 
we can adapt the number of sections to the period of the orbit. 

(6) We move the interval currentEpsEsinge and proceed to the next step of the 
continuation. 

The continuation starts with incrementSize = 10“® and the size (diameter) of 
the h-set Xq of order 10“® and both of these parameters vary throughout the proof. 
If any step of proveExistenceOf Orbit fails - for example Newton’s method does 
not converge or there is no covering between the h-sets, the algorithm will try to 
redo all the steps for a decreased incrementSize and proportionally decrease the 
size of the initial h-set. If the algorithm keeps succeeding, the program will try to 
increase incrementSize and the diameter to speed up the continuation procedure. 
The theorem is proved when bounds of currentEpsEange pass the bounds of e we 
intended to reach. Values of incrementSize for several different currentEpsEsinge 
can be found in Table along with periods of the periodic orbit and amounts of 
sections given by pm_count. 


CurrentEpsEange 

incrementSize 

period 

pm count 

[0.0014933550,0.001499146] 

5.7918161 X 10-® 

[201.35884, 207.17313] 

179 

[0.001,0.001001] 

Hr® 

[283.37351, 292.02862] 

212 

[4.9947443, 5.0200138] x lO^* 

2.5269501 X 10^® 

J521.07987, 557.55718] 

301 

[1.5057754,1.5132376] x 10-^ 

7.4621539 x lO^’’ 

J1593.3303, 1846.4787] 

671 


Table 2. Sample values from the validated continuation proof of 
Theorem |1.2[ As one can see, the period increases significantly 
as e —>■ 0, making it necessary to introduce more sections and 
lengthening the computations. 


4.3.1. Further continuation. We have decided to stop the validated continuation at 
e = 0.0015. Above that value our continuation algorithm encountered difficulties in 
its nonrigorous part, and needed many manual readjustments of the continuation 
parameters. As we later checked with MATCONT, this seemed not to have been 
caused by any bifurcation, so, most likely, it was just a defect of our ad-hoc method 
of continuing approximations of the periodic orbit by computation of Poincare maps 
between sections. Nonrigorous continuation methods implemented in continuation 
packages such as MATCONT are based on approximation of the orbit curve by 
Legendre polynomials and seem more reliable than our approach. Such a good 
nonrigorous approximation with a large number of collocation points would be 
enough to have a rigorous part of the continuation based on Poincare maps succeed, 
making further continuation only a matter of computation time. We did not 
implement it though, as we have decided that we are satisfied with how wide our e 
range is. By Theorem |1.3| we have already reached the value where the standard 
interval Newton-Moore method applied to a sequence of Poincare maps succeeds, 
and we think it is clear that a proof for higher values of e will pose no significant 
theoretical or computational challenges. 


4.4. Proof of Theorem 1.3 , Recall the interval Newton-Moore method for finding 


zeroes of a smooth map T : 
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Theorem 4.1 (The interval Newton-Moore method 2 p9pl] ). Let X = 6^], 

J- : K" — >■ M" be of class and let xq G X. Assume the interval enclosure of 
DJ-{X), denoted by [DJ-{X)\ is invertible. We denote by 

Af{xo, X) := -[DX{X)]-^X{xo) + Xq (130) 


the interval Newton operator. Then 

• if Af {xq, X) C intX, then the map T has a unique zero a;* G X. Moreover 
x* G Af{xo, X). 

• IfAf{xo, X)r\X = %, then X[x) 0 for all x G X. 

We applied the interval Newton-Moore method to a problem of the form ([^ given 
by the sequence of 179 Poincare maps obtained from the last step of the continuation 
procedure described in Subsection |4.3[ i.e. the step, where currentEpsRcinge 
contains 0.0015. Let i?max(0,7’) denote an open ball of radius r centered at 0 in 
maximum norm. We obtained the following inclusion 

M (o, (0,10-6)) C B„,ax (0,4.7926638 x lO'^^), (131) 

which, by Theorem |4.H implies the existence and local uniqueness of the periodic 
orbit. 


Remark 4.2. We report that we have succeeded with a verihed continuation based on 
the interval Newton-Moore method for the whole parameter range of Theorem |1.2[ 
that is e € [1.5 x 10“"^, 0.0015]. Although we got a little extra information on the 
local uniqueness of the solution of the problem ([^ , we have decided to discard this 
result, as it was vastly outperformed in terms of computation time by the method 
of covering relation^ It seems that the sequential covering process in the method 
of covering relations benefits more from the strong hyperbolicity than the interval 
Newton operator, hence allowing to make wider steps in the parameter range for 
such type of problems. However, for ranges of higher values of e the interval Newton- 
Moore method was only several times slower than the one of covering relations (e.g. 

in its 


7 times in the range [0.001,0.0015]), so we decided to state Theorem 1.3 


current form to show that we have achieved a parameter value where the more 
widespread tool is already adequate to the task. 

4.5. Technical data and computation times. All computations were performed 
on a laptop equipped with Intel Core i7 CPU, 1.80 GHz processor, 4GB RAM and 
a Linux operating system with gcc-5.2.0. We used the 568th Subversion revision of 
the CAPD libraries. The programs were not parallelized. 

Verification of assumptions of Theorem o took 236 seconds. Over 95% of the 
processor time was taken by verification of isolation for the chains of isolating 
segments. 

Proofs of Theorems |1.2| and |1.3| were executed by the same program. The validated 
continuation in Theorem |1.2| was the most time consuming part - it took 4153 


seconds. Theorem 1.3 is formulated for a single parameter value; the proof here was 
instantaneous - it finished within 2 seconds. 

We remark that the successful attempt to check the assumptions of Theorem o 
also for the range e G [10 1.5 x 10 ^] (119 seconds) saved us a lot of computation 

time. In theory we could have tried to use a validated continuation approach like in 
Theorem 1 1.2 for this range. We tried it later for a subrange e G [1.1 x 10-“^, 1.5 x 10“^] 


(for the whole range execution of Newton’s method for the problem ([^ within the 
nonrigorous part of the continuation algorithm failed due to enormous sizes of 
matrices to invert) - it took 2571 seconds, that is over 20 times longer. This 


^Substituting the interval Krawczyk operator for the interval Newton operator did not resolve 
this issue, i.e. did not allow for greater widths in the parameter steps. 
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indicates that construction of isolating segments around slow manifolds can be 
a valuable tool for proofs for “regular” parameter ranges (i.e. not including the 
singular perturbation parameter value) in systems with a very large separation of 
time scales. 


5. Concluding remarks 

We proved the existence of a periodic orbit in the FitzHugh-Nagumo equations 
for a range of the small parameter e S (0, 0.0015]. We also showed that the range is 
wide enough to reach its upper bound with standard validated continuation and 
rigorous methods. We hope that by further development of methods aimed at 
rigorous computations many classical results from singular perturbation theory, such 
as 

• a proof of existence of the homoclinic orbit in the FitzHugh-Nagumo system, 

• proofs of existence of periodic and connecting orbits for fast-slow systems 
with higher-dimensional slow manifolds, 

• proofs of uniqueness and stability of the waves, 
should be achievable in such explicit parameter ranges. 
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Note added in proof 


After the first preprint of this article was released, two other manuscripts on 
the existence of waves in explicit ranges of e appeared online. Matsue proved the 
existence of the periodic orbit, the homoclinic orbit and the heteroclinic loop for 
the FitzHugh-Nagumo system, all without further validated continuation 28 . The 


methods of Matsue are computer-assisted and similar to the ones described in this 
paper, though they seem to follow GSPT techniques more closely and put a bigger 
emphasis on tracking of invariant manifolds, at the cost of narrower ranges of e. 
The second work is a dissertation of AC 10 , which contains the results from this 


paper on the periodic orbit and extends them to give a proof of existence of the 
homoclinic, also without further continuation. 
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